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Ground states of interacting QFTs are non-gaussian states, i.e. their connected n-point correlation 
functions do not vanish for n > 2, in contrast to the free QFT case. We show that when the ground 
state of an interacting QFT evolves under a free massive QFT for a long time (a scenario that can 
be realised by a Quantum Quench), the connected correlation functions decay and all local physical 
observables equilibrate to values that are given by a gaussian density matrix that keeps memory only 
of the two-point initial correlation function. The argument hinges upon the fundamental physical 
principle of cluster decomposition, which is valid for the ground state of a general QFT. An analogous 
result was already known to be valid in the case of d = 1 spatial dimensions, where it is a special case 
of the so-called Generalised Gibbs Ensemble (GGE) hypothesis, and we now generalise it to higher 
dimensions. Moreover in the case of massless free evolution, despite the fact that the evolution may 
not lead to equilibration but unbounded increase of correlations with time instead, the GGE gives 
correctly the leading order asymptotic behaviour of correlation functions in the thermodynamic and 
large time limit. The demonstration is performed in the context of bosonic relativistic QFT, but 
the arguments apply more generally. 

PACS numbers: 03.70.+k, ll.10.-z, 05.30.-d, 05.30.Ch 


I. INTRODUCTION 

The nature of equilibration in quantum statistical physics [1-3] is still an open problem. Imagine that we prepare 
an isolated quantum system in the ground state of some hamiltonian and let it evolve unitarily under a different 
hamiltonian. This can be implemented by an instantaneous change of the hamiltonian of the system, a protocol 
known as Quantum Quench [4]. Since the initial state is a pure state, the system will always remain in a pure state 
and its evolution will be periodic or quasi-periodic, i.e. if we wait long enough it will return to the initial state or 
as close to it as desired. In the contrary, if we consider subsystems of the whole system, since they are not isolated, 
they are described by a reduced density matrix, which may correspond to a statistical ensemble. The period of 
quantum recurrences typically diverges with the system size, so that considering first the thermodynamic limit (i.e. 
infinite system size) and then the large time limit [5] it is possible that the system exhibits stationary behaviour, 
at the level of its subsystems and local physical observables. This behaviour however is not necessarily thermal; in 
fact one-dinrensional integrable systems exhibit equilibration to non thermal ensembles [6] . Progress in experimental 
techniques has made it possible to observe quantum non-equilibrium dynamics [7-15] and the occurrence of stationary 
behaviour different from thermal at large times [8, 16]. Theoretical attention has mainly focused on the equilibration 
mechanism in one-dimensional integrable systems [4, 6, 17-36] (see Ref. [1-3] for general reviews) [118]. 

In general (although not without exceptions), for non-integrable extended quantum systems it is expected that 
the stationary behaviour can be described by a thermal ensemble (Gibbs or Canonical ensemble) [37-42]. This 
expectation is based on the statistical principle of entropy maximisation under the constraint of energy conservation. 
In this scenario equilibration is equivalent to standard thermalisation. If we consider integrable systems, on the other 
hand, the picture changes dramatically: the number of local conserved charges is infinite and, due to these extra 
constraints, thermalisation is prohibited. The obvious generalisation is to construct an ensemble that contains all 
infinite conserved quantities [43], but still being economic (in the sense that not all projections into energy eigenstates 
are to be included, as this would amount to keeping all unnecessary information about the initial state). In the 
seminal paper [6] a Generalised Gibbs Ensemble (GGE) that includes all local conserved charges was proposed to 
describe the long time stationary behaviour of one-dimensional integrable systems. This conjecture has been proved 
correct in many cases, in most of which the evolution is performed under a free or effectively free theory, but also in 
special cases of quantum quenches in genuinely interacting integrable theories [5, 18, 21-23, 37, 44-67]. Even in cases 
where failure of the GGE in its original form was observed (some cases of interacting integrable models possessing 
bound-state eigenstates [68-73]), an alternative version of the GGE, sometimes called “Quench Action” ensemble [30] 
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or “overlap Thermodynamic Bether Ansatz” remains still valid [70, 72], while the insertion of quasi-local conserved 
charges [74-78] into the definition of the original GGE has been suggested as a suitable modification. 

It is tempting to consider the generalisation of the GGE idea in higher spatial dimensions. An analogous gen¬ 
eralisation of results for the transport properties of one-dimensional quantum systems to d > 1 spatial dimensions 
has been recently performed [79, 80]. Moreover new applications of the gauge-gravity duality to many-body physics 
[81-85] are also suggestive in this direction. Given that the GGE hypothesis refers to integrable models, extending its 
validity to higher dimensions would be restricted by the same restrictions imposed by integrability itself. According 
to the Coleman-Mandula theorem [86] and its extensions [87, 88], integrability imposes strong constraints to QFTs in 
dimensions d > 1 so that any integrable massive relativistic bosonic field theory in d > 1 must be free. Therefore we 
focus on the free case of the Klein-Gordon theory in d > 1, where the GGE conjecture is plausible despite the dimen¬ 
sionality. More generally, in the context of Quantum Field Theory, the introduction of the AdS/CFT conjecture [89] 
showed the first examples of seemingly integrable theories in spatial dimensions d > 1 [90] : Af = 4 Super-Yang-Mills 
(SYM) and Af = 6 Super-Chern-Simons theories [91], that are conformal and maximally supersymmetric, respectively 
in d = 3 and d = 2. The check of validity of the GGE hypothesis in such supersymmetric theories is left for future 
investigation. It should also be mentioned that the GGE may be useful in higher dimensions even for systems evolving 
under an interacting hamiltonian, not as an exact description of the large time limit (since there are no exact local 
conserved charges), but rather as an approximate description of pre-thermalisation, i.e. of the transient stage of the 
evolution that is expected to gradually lead to thermalisation [92, 93]. 

In the case of free evolution in one-dimensional bosonic models, validity of the GGE has been demonstrated in 
various examples of gaussian initial states (i.e. for free-to-free quantum quenches) or equivalently states with factorised 
charge product expectation values. In Ref.s [27, 34] it was shown that, neither gaussianity, nor factorisation of charge 
products in the initial state are necessary for the validity of the GGE hypothesis in d = 1. Instead, in Ref. [34] it was 
shown that the GGE conjecture is valid for a general interacting-to-free quantum quench in a one-dimensional bosonic 
system and the cluster decomposition principle, which is a fundamental physical requirement for the ground state of 
any physical hamiltonian, was identified as the sufficient and necessary property of the initial state that ensures the 
validity of the GGE. It is not a priori obvious whether dimensionality plays a crucial role in the above conclusion or 
not. In this paper we investigate this question. We find that in the case of massive free post-quench hamiltonian in 
d = 2 or 3, the GGE conjecture is still valid, i.e. the system equilibrates and the GGE describes correctly its stationary 
behaviour at large times. In the case of massless free post-quench hamiltonian on the other hand, the system does not 
always equilibrate but instead its field correlations keep increasing unbounded with time. In particular, in d = 2 the 
system never equilibrates, no matter what the initial state is, while in d = 3 it may or may not equilibrate, depending 
on whether the initial correlations decay exponentially or algebraically with a sufficiently low exponent, respectively. 
When equilibration does occur, the GGE turns out to be valid. We also observe that even in the cases where the 
system does not equilibrate, the leading asymptotic behaviour of long time correlations in the thermodynamic limit 
depends on only the two-point initial correlation function. All other information about the initial state contributes 
only subleading corrections in the large time and thermodynamic limit, although these corrections may in fact diverge 
with the system size or time (of course slower than the leading order expressions). The GGE predicts correctly only 
the leading order behaviour in this limit, so it is still applicable but in a weak sense. The economy property that the 
leading large time behaviour depends on only partial information about the initial state, which is characteristic of the 
GGE, seems therefore to be a more general feature of unitary evolution of extended quantum systems. 

The paper is organised as follows: In Section II we introduce the general concepts of equilibration in integrable 
models, the GGE hypothesis, the role of locality and of the double limit, thermodynamic limit and large time limit. 
Next we calculate the large time form of correlation functions after an interacting-to-free quantum quench and compare 
them with the GGE predictions. We first focus on the case of massive evolution in Section III, while in Section IV 
we study separately the case of massless evolution. Conclusions and outlook are discussed in the last Section V. 
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II. EQUILIBRATION AND GGE 
A. Definition of the model and quench protocol 

We consider the Klein-Gordon field theory in Minkowski space-time with spatial dimension d > 1 and defined in a 
large box of volume L d with periodic boundary conditions. This system is described by the hamiltonian 

H = ^ J d d x ^7 t 2 (x) — (V0( x)) 2 +m 2 (/) 2 (x)j 
= 2 'y . (fikK-k T 

k 

= ^ Eko\,a k + const, (1) 

fc 

where a],, a*, are the creation and annihilation operators and Ek = Vfc 2 + rn 2 is the relativistic dispersion relation. 
For completeness, we report standard definitions for this theory that are used throughout the subsequent calculations 
in Appendix A. 

We will eventually consider the thermodynamic limit L —> oo, since equilibration of the large time expectation 
values of local observables is only possible after this limit has been taken (the evolution of observables of any finite 
system exhibits quantum recurrences, i.e. periodic or quasi periodic behaviour, as will be discussed in more detail 
below). We will first focus on the massive case m ^ 0. This hamiltonian will play the role of the post-quench 
hamiltonian in our problem. 

The initial state |f2o) is assumed to be the ground state of another hamiltonian, which is left arbitrary for the 
moment, except that it is assumed to be symmetric under translations and rotations in space (and therefore the 
initial state is also translationally and rotationally symmetric). We will later distinguish two cases: when the initial 
state is the ground state of another free hamiltonian (with a different dispersion relation e.g. a different mass) or the 
ground state of an interacting hamiltonian. In the first case the initial state is gaussian, while in the second case it is 
not. In both cases the evolution is under a gaussian hamiltonian. 

Gaussian states are defined by their characteristic property that their only non-vanishing connected correlation 
functions are the one- and two-point functions. All other multi-point correlation functions in such states can be 
expressed diagrammatically in terms of the latter. More explicitly, choosing the one-point function to be zero (which 
can be done by a simple redefinition of the field cf> so as to subtract its expectation value on the state), the multi-point 
correlation functions can be written in terms of the two-point function by means of Wick’s theorem , i.e. by summing 
over all possible pair contractions. Examples of gaussian states are the ground states of free (quadratic) hamiltonians, 
but also their coherent and squeezed coherent states. Non-pure (mixed) gaussian states include thermal ensembles 
of free hamiltonians, as well as the GGE defined below, since all these can be expressed as exponentials of operators 
quadratic in the field </> and its canonically conjugate field 7r. It should be emphasised that gaussianity refers to the 
particular choice of fields whose correlation functions we are interested in. In the present problem, ‘gaussian initial 
state’ means that the above property holds for the correlation functions of the local fields cf> and 7r which are also the 
one in terms of which the post-quench hamiltonian H is free. 

We will also assume that the initial state is symmetric under the transformation (f> —» —0. This means that the 
initial one-point function vanishes ( <fi(x))o = (STo |<^(aj) |A2o) = 0. In the opposite case, which is valid for example 
when the pre-quench ground state exhibits spontaneous symmetry breaking, the one-point function develops a trivial 
oscillatory behaviour (Appendix B). 


B. The GGE conjecture, locality and extensivity aspects 


According to the GGE conjecture, the expectation value of any local observable (or subsystem of the system) in 
the thermodynamic and large time limit is given by a statistical ensemble, described by the GGE density matrix 


PGGE = exp 



( 2 ) 


where Qi are the local conserved charges of the integrable system under evolution and each /3j is the Lagrange 
multiplier associated to Qi accounting for the respective conservation law. The values of /3j are determined by the 
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self-consistency condition that the values of the charges in the GGE ( Qi)gge (being local observables themselves) 
are the same as their values in the initial state (Qi)o 

(Qi)GGE = T ffi PgC *> = (Ool W) = (Qi) 0- (3) 

Trj Pgge} 

This is the only information about the initial state that is contained in the GGE. 

The GGE in its original form is supposed to contain all local conserved charges Qi and only those [5, 6]; no other 
constant of motion is included [119]. However it is typically written in terms of occupation number operators of 
all independent modes in which the integrable system can be decomposed, which are non-local operators. In Klein- 
Gordon theory these are the occupation number operators of each of the momentum modes n(k) = a^ak which are 
conserved operators and manifestly non-local. In terms of these operators the (non-normalised) GGE density matrix 
is 


Pgge = exp 
with ( n(k)) GGE 


( k )n(k) 

k 

= ( n(k)) 0 ■ 


( 4 ) 


In free field theories, like in the present case, this is an equivalent alternative form because, the set of independent 
local charges is a linear combination of the occupation number operators, provided there are no convergence problems 
arising in the continuum limit [78]. In fact in most demonstrations of validity of the GGE conjecture, the latter 
is written in the form (4). The locality requirement is related to the expectation that in a statistical ensemble 
the (generalised) “internal energy” Egge should be an extensive quantity, in order to ensure the extensivity of 
thermodynamic quantities of macroscopic subsystems. 

However in some cases of genuinely interacting integrable models, the GGE in its original form in terms of local 
charges (2) has been proven to fail to describe the large time limit [68-73], while in the form (4) and taking into 
account all independent modes it is still correct. Moreover (4) gives correct predictions even in cases where it is not 
possible to express it in terms of local charges and even if there do not exist local charges at all, like in the case of 
condensed matter models confined in external potential traps or models with long-range harmonic couplings [94, 95]. 
In such cases the GGE internal energy may be non-extensive. 

In the present problem the GGE internal energy corresponding to (4) is 

Egge = - log Zgge + Sgge = f3{k)(n(k))GGE = — ^ dfog/g^ ’ ( 5 ) 

where Sgge is the GGE von Neumann entropy 

Sgge = — ^{p'gge ^°&p'gge} ( 6 ) 

with 


the normalised GGE density matrix and 

Zgge = Tr pgge-, (8) 


Pgge = 


Pgge 
Zgge ’ 


( 7 ) 


the GGE partition function. Expressing n(k) in terms of <f>(x, t ) and n(x, t) = <fi( x , t) fields and using the translational 
invariance of the initial state, we finally find that 


Egge = ^ J d d r [ ZMr )^ 0 , 0 ) ( 0 , r ) + iD 2 (r){C^’ 1] (0, r) - C ^' 0 ) ( 0 , r )) + D 3 (r)C i Q 1 ’ 1) (0,r)] 


(9) 


c'o°’ 0) ( a: , y) = (<!>(x)<l>(y))o, 

(10) 

C'o 0 , 1 \x,y) = {(/)(x)TT(y))o, 

(11) 

Co 1 , 0 ) (x,y) = {Tr(x)<j)(y)) 0 , 

(12) 

C ( o^\x,y) = (n(x)ir(y))o, 

(13) 


where 
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are the initial correlation functions and 


£>i(r) = ^ 

[ d d k e ik ' r j3(k)Ek, 

(14) 

D 2 (r) = ^ 

j d d k e ik ' r p(k), 

(15) 

D 3 (r) = ^ 

[ d d k e ik - r p{k)/E k . 

(16) 


Therefore, if f3(k),E k and ( n(k))GGE are analytic functions in a neighbourhood of k = 0 (and decay sufficiently fast 
at large k), then Di[r), i = 1,2,3, all decay exponentially fast with the distance r and the internal energy Eqge is 
extensive. Note that in the massless case m = 0 the dispersion relation E k = |fe| is not analytic at k = 0, therefore 
the internal energy is a non-extensive quantity in this case. 

In any free field theory the GGE density matrix in the form (4) is a gaussian ensemble with respect to the bosonic 
field 0, i.e. the multi-point correlation functions of <f> satisfy Wick’s theorem [120]. This is because it is a function of 
the momentum occupation number operators n(fc) which are quadratic in terms of the bosonic field and its canonical 
momentum. Therefore if the GGE conjecture is valid, it equivalently means that the gaussian evolution erases the 
memory of non-gaussian initial held correlations, i.e. that connected correlation functions decay with time in the 
thermodynamic limit. 


C. Validity of the GGE and the cluster decomposition argument 

We now outline our method to check the validity of the GGE conjecture based on the cluster decomposition principle. 
We have to calculate the evolution of local observables and, provided they equilibrate, to compare their large time 
stationary values with the GGE predictions. We will focus on the correlation functions of the bosonic held cf>(x : t). 
For the two-point function, testing the GGE is trivial, since, provided it equilibrates, its stationary value is given in 
terms of the conserved occupation numbers of the momentum modes n(k), which by dehnition of the GGE are exactly 
the information contained in it. The hrst nontrivial test comes from the study of the four-point function, since its 
stationary value depends also on the initial correlations of n(fc), i.e. on (n(k)n(p ))o which are also conserved by the 
evolution but their values are not taken into account in the construction of the GGE. The GGE prediction for their 
values is simply the factorised value {n(k))o{n(p))o, since it is a gaussian ensemble and Wick’s theorem is valid in it. 
In the quantum quench case, if the pre-quench hamiltonian is also free, then Wick’s theorem is valid in the initial state 
too and therefore ( n(k)n(p)) 0 = {n(k))o(n(p))o and the GGE prediction is correct. If the pre-quench hamiltonian is 
not free, then (n(k)n(p ))o does not factorise, however in the previously studied one-dimensional post-quench-massive 
case [34] the GGE is still valid. The reason is that the contribution of connected correlations 

(n(fe)n(p))o - (n(k)) 0 {n(p)) 0 (17) 

to the large time value of the four-point function turns out to be a finite size correction, which vanishes in the 
thermodynamic limit L —> oo. To see this we express the large time limit of correlation functions in terms of initial 
correlation functions and realise that they enter as a weighted spatial average over the whole system. The initial state 
after a quantum quench, being the ground state of some physical hamiltonian, satisfies the cluster decomposition 
principle, meaning that at large separations all of its multi-point correlation functions factorise to the fully contracted 
expression. As a consequence, when averaging them over the whole system, the only surviving information is the 
initial two-point correlation function, or equivalently the quantities (n(fe))o, which information is precisely contained 
in the GGE. This argument works similarly for the large time limit of any other multi-point correlation function of 
order higher than four. The above reasoning is diagrammatically represented in the diagram of Fig. 1. 


D. Finite size, recurrences and the order of large time and thermodynamic limit 

The order of the two limits t —»• oo and L —oo in the above method should be clarified. At finite volume the system 
exhibits recurrences (also called revivals), i.e. all correlation functions are periodic or quasi-periodic functions of time 
with a recurrence period that increases with L [2]. This is because the energy spectrum of a finite system is discrete 
with the spacing between two successive energy levels being of order 1/L, due to the quantisation of momenta. If 
the energy eigenvalues are commensurate, i.e. they are all integer multiples of a fundamental frequency, then the 
evolution is exactly periodic. This is the case in Id massless relativistic systems and the recurrence period is equal to 
the system size L , since the dispersion relation is E = |fe| = 27r|n|/L, with n integer. If instead the energy eigenvalues 
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FIG. 1: Diagrammatic representation of the method. 


are incommensurate, then revivals are partial and the evolution is quasi-periodic, meaning that the system will return 
as close to the initial state as we want, if we wait long enough for a sufficiently approximate revival. 

From the above we see that the limit t —» oo does not exist at finite L. Instead, if the thermodynamic limit L —> oo 
is taken first, then the evolution is no longer periodic (the recurrence period diverges) and it is now possible for the 
system and its correlation functions to exhibit stationary behaviour as t —> oo. This means that the two limits do not 
commute 


lim lim ^ lim lim , 

t—¥ oo L—t OO L—t OO t—¥ oo 


(18) 


as one may exist while the other does not. However let us consider the “long time average” values of correlation 
functions at finite volume L, defined as 


C (n \x i, ...,x n ) = lim T 

T —too 


dtCi n \x 


(19) 


where c[ n ' > is the correlation function at time t (and at finite volume). In fact if the evolution is strictly periodic 
(not quasi-periodic), then it is sufficient to average over a single period, which as mentioned above is an increasing 
function of L. Obviously the long time average values depend explicitly on L. If the system becomes stationary in the 
limit Hindoo lirn^oo, then the stationary values of correlation functions are identical to the thermodynamic limit of 
the long time averaged values, i.e. 

r T 

lim lim T -1 / dt c[ n \x\, ...,x n ) = lim lim C^ n \xx n ) (if equilibrated). (20) 

L—tooT—too J t—too L—too 

Therefore the two limits can be effectively permuted if we employ time averaging. This has the technical advantage 
that long time averages in finite systems can be calculated easier. 
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and its long time average value (dashed lines) for system sizes L that increase successively by a factor of two (darker grey 
lines correspond to larger sizes L), along with their asymptotic behaviour in the thermodynamic limit L —> oo (red full and 
dashed lines). Top: the case of equilibration. Bottom: the case of unbounded increase. Notice that in the first case the time 
averages converge as L —> oo to the stationary value, while in the second case the unbounded increase with the time t in the 
thermodynamic limit is reflected in the finite size system as unbounded scaling of the time averages with increasing L. 


If, instead of strict equilibration, the system exhibits bounded persistent oscillations in the limit lirn^oo lim^^oo, 
these will obviously be lost after time averaging. If on the other hand there is neither equilibration nor equilibration 
on average, but instead correlations keep increasing indefinitely for large times, then the long time average will be 
divergent with L in the thermodynamic limit (Fig. 2). In general, equilibration means that the evolution of observables 
in finite size systems exhibits a “plateau” (typically accompanied by rapid finite size oscillations) during a long time 
window between two recurrences, that becomes infinite in the thermodynamic limit [2, 3]. Instead, if the system does 
not exhibit such a plateau, then the time averaged observables may increase, when the system size (and therefore 
the recurrence period) increases, and finally diverge in the thermodynamic limit. We can therefore study the scaling 
behaviour of the time averaged correlation functions as functions of L, in order to determine the scaling behaviour of 
the correlation functions as functions of time in the thermodynamic limit and determine whether they equilibrate or 
not. 


III. MASSIVE NON-INTERACTING POST-QUENCH HAMILTONIAN 

Having outlined our method and clarified certain important points, we are now ready to check the validity of the 
GGE, first focusing on the case of massive post-quench hamiltonian. We will therefore compare the long time averages 
of correlation functions of the held 4>(x,t ) with the GGE predictions. The calculation is presented in full detail in 
the Appendices. For the two-point correlation function C\ (x\,X 2 ) = (<j){xi\t)(j>{x 2 \t)) we find that its long time 
average value is (Appendix C) 

C^\x 1 ,x 2 ) = & 2 \x i - x 2 ) 

e ifc .(* 1 -* 2 )( 1 + 2 ( n ( fc )) 0 ). ( 21 ) 

i.e. it is expressed in terms of the momentum distribution ( n(k))o in the initial state, which is conserved under the 
evolution. The above expression can be written in terms of coordinate space initial correlation functions as 

C (2) ( r ) = * (<7< 2 > ( r ) + J d d s H(s)C < 2) (r - s)^j , ( 22 ) 


d d k 1 
(2^p2 E t 


where we defined 


f d d k 




H(r) = 


(2n) d El ’ 


(23) 






















and 


C'o 2) ( r ) = ( 7r ( r i 0)7r(0,0)) 

d d 


dt\ dt 2 




(24) 


t\ —^ 2—0 


This should be compared to the GGE prediction, which is (Appendix F) 


c gge ( x i - x z) 


d d k 


( 2 n) d 2 E k 


e ik (x 1 -x 2 )^ _|_ 2(n( k )) GGE ). 


(25) 


i.e. it is expressed in terms of the momentum distribution in the GGE ( n(k))QQE . The two expressions (21) and (25) 
are trivially identical, i.e. 


due to the defining property of the GGE that 


0 (2 Hr) = C^%(r) 


(n(k)) G GE = ( n(k )) 0 • 


(26) 


(27) 


Next, for the four-point correlation function C^ 4 \x 1 , x 2 , x 3 , x 4 ) = ((/>(x 1 ;t)<j>(x 2 ; t)(j>{x 3 \ t)<p(x 4 -, t)), the long time 
average value is (Appendix D) 


C (4 \x 1 ,x 2 ,x 3 ,x 4 ) 


1 

2 


d d k d d p 1 
(27T) 2d 4 E k E v 


S(k,p; x 1 ,x 2 ,x 3 ,x 4 ) 


X 


[n(k)n(p)) 0 + (n(k)) 0 + - 


(28) 


where 


5(fc,p;® 1} ® 2 ,*3,*4)= e lk < X2 - Xl)+tp ^- X3) . 

all perm.s 
of 1,2,3,4 


In terms of initial correlation functions in coordinate space, the above can be written as 


(29) 


C ( 4 ) ( a ; i ,® 2 ,® 3 ,® 4 ) = 


1 

32 


E 

all perm.s 
of 1,2,3,4 


L d 


d d sd d r 

L 2d 


Cq 4) (s, s + Xi 


X 2 ,r,r + x 3 - ® 4 ) 


+ J d d s H(s')Cq 4 \s,s' + s + ®i - x 2 ,r,r + x 3 - ® 4 ) 

L d 


+ J d d r' H(r'){Co‘\s, s + ®i - x 2 ,r,r' + r + x 3 - ® 4 ))* 

L d 


+ JJ d d s'd d r' H(s')H(r') Cq^{s, s' + s + x 3 

L d 


x 2 ,r,r' + r + x 3 - ® 4 ) , 


(30) 


where 


Cq 4 \x 1 ,x 2 ,x 3 ,x 4 ) = (tt(x 1 )tt(x 2 )(I)(x 3 )(I)(x 4 )) 0 (31) 

and 

C < 0 4 \x 1 ,x 2 ,x 3 ,x 4 ) = (ir(xi)TT(x 2 )n(x 3 )Tr(x 4 )) 0 ( 32 ) 
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while the function H(x) has been defined in (23). The last relation shows that the long time average of the four-point 
function C ^ can be expressed as a weighted spatial average of initial correlation functions C ^, C'q 4 '* and C g 4 '* ■ 
Notice that (28) depends not only on the initial momentum occupation numbers, but also on their initial correlations 
( n(k)n(p))o. In contrast, the GGE prediction is given in terms of solely the momentum occupation numbers ( n(k)) 0 . 
Indeed, since the GGE for a free model is a gaussian ensemble, Wick’s theorem applies and therefore the correlations 
of the occupation numbers factorise 


(n(k)n(p)) gge = ( n(k)) G GE(n{p)) G GE , 
thus leading to the prediction (Appendix F) 


C GGe( X 1, x 2,X 3 ,X 4 ) 

1 f d d k d d p 1 


2 J (2t r) 2 4 E k E p 


S{k,p;x 1 ,x 2 ,x 3 ,x 4 ) 


x ( i n {k)) gge{ n{p ))gge + {n{k)) GG E + ^ 


(33) 


(34) 


or simpler, in coordinate space 


Cgge( x i> x 2, x 3,x 4 ) = Cggei^I’ x 2)Cgge( x 3> x i) 

+ [2h3] + [2h4], (35) 

A. Non-interacting pre-quench hamiltonian 

From the comparison of (28) and (34), an obvious case for which the GGE prediction is valid, is when the correlations 
of charges factorise also in the initial state 


(n(fc)n(p)) o = (n(k))o(n(p)} 0 , (36) 

i.e. when Wick’s theorem is valid in the initial state, as is the case of a ground state of a pre-quench hamiltonian that 
is also free (i.e. quadratic) in terms of the bosonic field (j), though with different physical parameters (e.g. different 
mass). Therefore for a free-to-free quantum quench the GGE conjecture is always valid. This is because the initial 
state is then gaussian, like the GGE itself, and therefore the consistency of the two-point function expressions (21) and 
(25) is sufficient to ensure consistency for the four-point function and, following the same reasoning, for all multi-point 
correlation functions. 


B. Interacting pre-quench hamiltonian 

We now turn to the more interesting case of non-gaussian initial states, such as the ground states of interacting 
hamiltonians. We will show that, even though the above factorisation does not hold, the spatial averages of the initial 
four-point correlation function entering in (30) are equal to the product of two-point correlation functions so that the 
GGE prediction is still correct. 

Let us focus on the first term in the integral of (30), which is essentially a spatial average of the initial four-point 
correlation function Cq 4 '* with respect to two coordinate variables s and r. From the general cluster expansion for 
the Cq 4) correlation function, taking into account the translational invariance of the initial state and our assumption 
that {(f{x)) 3 = (j) o, we have 


Cq 4) (s, s + aq - aq, r, r + x 3 - x 4 ) 

= C^ 2 ) (x i - a: 2 )G^ 2) (a;3 - x 4 ) 

+ Cq 2 \x - s)C 3 2 \r - s + x 3 - x 4 + x 2 - aq) 
+ C (2 \r - s + x 3 - x 4 )C 3 2 \r - s + x 2 - aq) 
+ CqcOT a:i — x 2l r — s,r — s + x 3 — aq) 


( 37 ) 
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We then see that substituting into the integral, only the first term in the above expansion, which is independent of s 
and r , survives. Indeed, due to the cluster decomposition principle all other terms give contributions that vanish in 
the thermodynamic limit. We therefore have 


lim 

L—too 



L d 


d d sd d r ^( 4 ) 
-^2d- C o 


( s , S + Xi 


= Cq 2) (Xl - x 2 )c{j 2 ) (x 3 - x 4 ). 


x 2 ,r,r + x 3 - x 4 ) 


(38) 


As shown in more detail in Appendix E, similar results hold for all other terms in the integral of (30). For the last 
term, for example, we find 


lim JJ ^^2 d r JJ d d s'd d r' H(s')H(r')C^\s,s' + s + xi - x 2 ,r,r' + r + x 3 - x 4 ) 

L d L d 

= JJ d d s'd d r' H(s')H(r') Cq 2 \x 3 - x 2 + s')Cq 2 \x 3 - x 4 + r'), 


(39) 


where the last expression is always finite as L —► oo for a massive post-quench hamiltonian, since the function H(r) 
decays exponentially with the distance r. 

Summing up all surviving terms, (30) reduces to 

C (4) ( x 4 ,x 2 ,x 3 ,x 4 ) = 

= ( c o ) (x 1 -x 2 )+ J d d s'H{s , )C^\x 1 -x 2 +s') 

all perm.s ' u 

of 1,2,3,4 

x ^Cq 2 ^ (x 3 — x 4 ) + J d d r'H (r , )C'o 2) ( x 3 - x 4 + r')j , (40) 


which, using our earlier result ( 22 ) for the two-point correlation function, can be identified with 

C ( - 4 \x ll x 2 , x 3l x 4 ) 

= i V] Cd 2 \x 4 — x 2 )C ( ' 2 \x 3 — x 4 ) 

o — 

all perm.s 
of 1,2,3,4 

= <7^(a;i - x 2 )C^ 2 \x 3 - x 4 ) + [2 q 3] + [2 q 4], (41) 

In this last form it is manifest that Wick’s theorem is valid for the long time average of the four-point function. 
Therefore comparing with the GGE expression (35) and since from (26) the GGE prediction for the two-point function 
is correct, we conclude that the GGE is correct also for the four-point function 

C (4) {x ll x 2 ,x 3 ,x 4 ) = C^l E (x 1 ,x 2 ,x 3 ,x 4 ), (42) 

always assuming that the post-quench hamiltonian is massive. 

Overall this means that, even when the initial state is a non-gaussian state, the contribution of connected initial 
correlations of the occupation number operators turns out to vanish in the thermodynamic limit, as a consequence of 
the cluster decomposition principle. Therefore as anticipated, the GGE conjecture is still valid in this case. 


C. Higher order correlation functions 

Higher order multi-point correlation functions can be calculated similarly and the above conclusions about equili¬ 
bration and validity of the GGE prediction hold the same. Indeed the stationary or long time average value of the 
2 n-point function C ( - 2n - ) ({*«}) can be expressed in terms of the expectation values of higher products of the momentum 
occupation number operators 


U(n{kj) + n {-kj) + 1 ) 

\i=i / 


( 43 ) 




11 


which in turn can be written as functions of the initial correlation functions 


n 

\j =i 



(44) 


where aj = 0,1 and we denote vf ' 1 = (j)k j , q = 4>kj ■ The result can be written in coordinate space as a 
weighted spatial average of the correlation functions 


+x 2j -i)4> {(7j \s j +x 2j ) \ 

\j=1 / r 


(45) 


which in the thermodynamic limit L —oo, due to the cluster decomposition property of the initial state, tends to 
the fully disconnected expression 


Y[(4>{x 2 j-i)4>(x2j))o, (46) 

l=i 

that is the gaussian result based only on the two-point correlation functions, precisely equal to the GGE prediction. 


IV. MASSLESS NON-INTERACTING POST-QUENCH HAMILTONIAN 

We now discuss the case of massless free evolution. This case should be studied separately, because the function 
H(r) does not decay exponentially with the distance as before, and this property is important in the calculation of 
the L-scaling of spatial averages of initial correlation functions (Appendix E). Instead of decaying exponentially, H(r) 
exhibits infrared divergences in d < 2, meaning that it diverges with L, while in d > 2 it decays as a power of the 
distance 


H(r) 


Id-2 


for \r\ oo and d > 2. 


(47) 


More specifically, its precise functional form in dimensions d = 1,2 and 3 is the following 

f (T/tt) - (r/2) if d = 1, 

H(r ) = H(r) = < log(L/r)/(27r) if d = 2, 

I l/(47rr) if d = 3. 


(48) 


This has as a direct consequence that the scaling behaviour of C^ 2 \xi, x 2 ) and C^{x\, x 2 , x 3 , X 4 ) is also different, 
depending on the dimensionality. Below we discuss separately the cases of d = 1,2 and 3. 


A. One dimension 

In d = 1, neither the two-point function nor the four-point function of the bosonic field become stationary in the 
large time limit, independently of the behaviour of the initial correlation functions Cq 2 \ Cq 2 \ Cq 4 \ Cq 4 ' 1 and C q 4 \ 
since both C^ 2 \xi : x 2 ) and C^\x i, x- 2 , X 3 , X 4 ) diverge in the thermodynamic limit L -A 00 . This observation is in 
agreement with what is known from [4] for initial states that belong to the class of conformally invariant boundary 
states (which are gaussian), where it is shown that the two-point function C\ (x\ 1 X 2 ) increases linearly with time. 
The divergence of & 2 ^(x\,x 2 ) with L, following the arguments in Section IID, reflects exactly this fact. However 
the d = 1 case is special in the sense that local physical observables in a massless one-dimensional held theory, i.e. 
in a Conformal Field Theory, are described by correlations of vertex operators instead of the free bosonic held <f> and 
those correlation functions do equilibrate (to values given by an effective thermal ensemble [4] in some special cases 
of initial states, or to the GGE in more general cases [45, 67]). 
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B. Two dimensions 


In dimensions d > 2 the ^-field correlations are physical observables themselves. In d = 2, due to the logarithmic 
divergence of H(r) with L, long time averaged correlation functions also diverge with L. Therefore, according to 
the above, there is no equilibration but instead correlations increase with time. Assuming that the initial correlation 
functions C g and Cq do not diverge with L themselves, which is always true for ground states of non-critical 
pre-quench hamiltonians, the leading order values of long time averages are dominated by the contribution of terms 
with the most H -weighted spatial integrals, i.e. 


and 


C^ir)^ 1 - J d d sH(s)C^(r-s), 

L d 

C ( 4 \x 1 ,x 2 ,x 3 ,x 4 ) 

= 32 ^ 2 5 E ff d d sd d rd d s'd d r'H(s')H{r') 

all perm.s , 
of 1,2,3,4 L d 

x C q 4 ) (s, s' + s + £Ci - x 2 l r 1 r' + r + x 3 - £c 4 ), 


(49) 


(50) 


which, using the same arguments based on the cluster decomposition principle (Appendix E), reduces to 


C (4) ( Xl i X 2 ? 3^3 j *£ 4 ) 

i - = ^° -L E j d d s'H{s')C^\x 1 -x 2 + s') 

all norm c ^ 


all perm.s , 
of 1,2,3,4 L d 


x J d d r'H{r')Co‘\x 3 — X 4 + r'). 

L d 


(51) 


( 2 ) 

If, on the other hand, Cq ' also diverges with L , which may happen when the pre-quench hamiltonian is critical too, 

and if it diverges also logarithmically, while the other initial correlation functions do not diverge with L, then the 

( 2 ) 

contribution of the terms involving Cq ' may be of equal importance with the rest, in which case we obtain expressions 
similar to (22) and (40), i.e. those holding in the massive-to-massive case. Similar observations of divergence of the 
two-point function in the special case of free-to-free quantum quenches in dimensions d = 2 have been made in 
Ref. [96], where it was shown that a mass quench in the Klein-Gordon theory leads to a two-point function that 
increases logarithmically with time. 


C. Three dimensions 

The d = 3 case is different from the above, because H[r ) decays with the distance as a power and therefore the 
answer to the question whether the correlation functions equilibrate or not depends also on the large distance decay of 
the initial correlations and in particular whether they decay sufficiently fast so that all spatial integrals are convergent. 
For the two-point function, if Cq ; decays exponentially, as is the case of a non-critical pre-quench hamiltonian, then 
equilibrates. The same is true if Cq 2 '* decays as a power but sufficiently fast, more specifically faster than r -2 , 
while if it decays as r ~ 2 or slower then c[ 2 ^ does not equilibrate but increases with t (i.e. C^ diverges with L ). 
Similarly for the four-point function Cj; 4 \ the large distance decay of Cq 2 ' 1 and Cq 4 ) determines if C) 4 '* equilibrates or 
not. Notice that, especially in the critical pre-quench case, it is possible that some correlation functions equilibrate 
while others do not. 


D. General remarks for d > 2 


According to the above, the scaling of long time averages of correlation functions as L —> 00 in dimensions d > 2 
(and therefore also the scaling of correlations as t 00 ) depends highly on the large distance behaviour of the initial 
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correlation functions Cg 2 \ Cq 2 \ Cq 4 ^, Cg 4 ^ and Cg 4 ^- This is determined by the Renormalisation Group (RG) analysis 
of the pre-quench hamiltonian. Note that the large distance behaviour of initial 7r-field correlations is related to that 
of the ^-field correlations, since n = <j>. In the case of the two-point function, for example, 7r-field correlations can be 
expressed, using general properties like space-time symmetries before the quench and Lorentz invariance, as a spatial 
derivative of the 0-field correlations 


M = -^cf(r). 


(52) 


Also note that the large distance decay of initial connected correlation functions of different orders is not independent 
of one another. If the four-point function C^ c decays as a power, the same must hold also for the two-point function 

Cq . Indeed, one contribution to the diagrammatic expansion of the two-point function involves the four-point 
function with two legs contracted forming a loop. This diagram decays as a power of the distance of the two 
uncontracted points and since the large distance decay is determined by the slowest decaying contribution, the full 
two-point function must decay as a power (Fig. 3). It should be mentioned that this argument is not perturbative. In 


C 4 = 


c,= 



FIG. 3: Pictorial representation of the argument that, if the connected four-point function decays as a power when any of the 
coordinates goes to infinity, then the same holds also for the two-point function. 


general we can distinguish two typical types of behaviour: initial connected correlation functions decay exponentially 
with the distance (as for non-critical pre-quench hamiltonians) or they decay algebraically (as for critical pre-quench 
hamiltonians). In the critical case it is also possible that correlation functions diverge with L themselves. 

In all of the above studied cases, both in d = 2 and d = 3, whenever all correlation functions equilibrate, their 
stationary values satisfy Wick’s theorem and are given by the GGE. Moreover, in the non-critical pre-quench case, 
even when correlation functions do not equilibrate, they are still expressed in terms of solely the two-point initial 
correlation functions (Cq and Cq ) at leading order in the large time and thermodynamic limit. Accordingly, the 
GGE predicts correctly (only) the leading order of long time averaged correlation functions in this case. Higher order 
correlation functions (Cg 4 c \ Cg 4 ! and c{, 4 ') give contributions next-to-leading order at large L, which may still diverge 
with L but slower. Therefore, even though Wick’s theorem (41) is not valid in this case, it is still applicable in the 
following weak sense: 


_ C (4) (a; 1 ,x?,x 3 ,x i ) _ z,->oo 

C( 2 )(:ei, x 2 )C( 2 \x 3: aj 4 ) + [2o3] + [2q4] 

Note that the difference between numerator and denominator may be divergent with L. Following the arguments 
of Section III C and Appendix E, this weak version of Wick’s theorem is expected to apply similarly to all higher 
order correlation functions. The above result is always true when the pre-quench hamiltonian is non-critical, but it is 
also true in the case of critical pre-quench hamiltonians for a wide range of power law exponents for the correlation 
functions. Even though a general classification of all possibilities is beyond the objectives of the present work, possible 
exceptions may arise only in marginal cases when the exponents of the algebraic decay of different correlation functions 
have special values, such that the contributions of connected four-point correlations functions survive at leading order 
in the large time and thermodynamic limit. 

It should also be emphasised that for any specific model, as the dimensionality d increases the deviations from the 
GGE become smaller and smaller until they finally vanish above some critical dimension. This is due to the scaling 
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of the convolutions involved in the expressions for the long time averaged correlation functions: both the function 
H(r) ~ l/r d ~ 2 and the initial correlation functions decay faster as d increases (even if they decay as a power law), 
so that the expressions above eventually become convergent. This can be seen easier in momentum space where the 
L —> oo divergences of the convolutions are determined by the order of the pole at k = 0. While the poles of H(k) 
and of the correlation functions are fixed, the phase space integration measure d d k ~ k d ~ 1 dk changes with d and 
eventually cancels the pole at k = 0. Therefore equilibration occurs and the GGE conjecture becomes correct for any 
interacting-to-free quench above some critical dimension. This is obviously the same argument as in mean field theory, 
which (both in the dynamical and in the equilibrium version) becomes exact above the critical dimension. 


V. CONCLUSIONS & OUTLOOK 

In this work we have studied the equilibration properties and validity of the GGE conjecture in higher dimensional 
bosonic QFTs, undergoing a quantum quench from a free or interacting hamiltonian to a free hamiltonian. If the 
post-quench mass is nonzero the system equilibrates, its correlation functions become gaussian i.e. satisfy Wick’s 
theorem and the GGE predicts correctly their stationary values. We can interpret this effect as loss of memory of 
initial non-gaussian correlations, due to evolution under a gaussian massive hamiltonian. From the technical point of 
view, considering the combination of thermodynamic and large time limit, it turns out that only partial information 
about the initial state is required to describe local observables and this information is precisely the initial two-point 
function averaged over space with a weight function H(r) dependent on the post-quench hamiltonian. If the post- 
quench hamiltonian is massless, despite the fact that the evolution may lead to unbounded increase of correlation 
functions, the GGE conjecture still gives correct predictions for the leading order correlations in the thermodynamic 
and large time limit, at least (but not only) in the case of non-critical pre-quench hamiltonians. Contributions of 
non-gaussian initial correlations may still increase with time or system size but slower, so that their strength relative 
to the leading order decays to zero. Moreover, above some model-dependent critical dimension, corrections to the 
GGE predictions vanish and the GGE conjecture becomes exact. Notice that our arguments are not sensitive to 
the particular quench under consideration, but refers to a large class of quenches and is expected to be valid for 
non-relativistic models as well. 

It is tempting to attempt to generalise our approach to the case of interacting evolution in order to test the hypothesis 
that in dimensions d > la suitable GGE, built with approximately conserved charges, describes the pre-thermalisation 
regime. This corresponds to the intermediate stage of the time evolution, after dephasing of quasi-particles takes place 
and before thermalisation due to their interaction is finally reached [29, 42, 110-117]. Our method relies on two key 
features of free evolution that are not generally valid in the interacting case: the solution of the equations of motion for 
the field operators and the locality of these fields, on which the cluster decomposition is based. However, perturbation 
theory suggests that these properties may be unaffected by the introduction of weak interactions, as long as: the 
dimensionality is d > 1, we are away from criticality and excitations have a finite maximum group velocity. These 
conditions are expected to guarantee the light-cone form of dynamics and the possibility to express the correlation 
functions in terms of local initial fields, which are sufficient conditions for the application of our method. In this 
perturbative approach, the approximate conserved charges would not be quadratic operators in terms of the local 
fields, which means that they would not be necessarily local either. We hope to come back to this problem in a future 
publication. 

Interesting questions arise also in the context of supersymmetric quantum field theories inZ) = d+ l>2 which are 
higher dimensional yet integrable. The natural candidates to explore the problem of equilibration in such theories are 
the AT = 4 Super-Yang-Mills (SYM) and Af = 6 Super-Chern-Simons theories (ABJM theory), which are conformal 
and in some sense also integrable [90] , though integrability in this case means partial knowledge of the spectrum of the 
theories, but not existence of infinite number of local conserved charges, as in one dimension. Indeed an infinite tower 
of non-local conserved charges appears in the scattering amplitudes of gluons (Yangian Symmetry [97]), but the role 
of these charges in determining the spectrum is not clear. There are two alternative approaches to study the quench 
dynamics in such models. One is based on a mapping of the quench problem to a boundary problem in a Euclidean 
slab geometry [4]. This was initially applied to Id Conformal Field Theory with a conformally invariant initial state 
which turns out to lead to effective thermalisation, while it was later generalised to perturbations of the latter states 
[45, 67] in which case a suitable GGE was shown to be correct. A natural extension to supersymmetric models would 
be to consider supersymmetric initial states [98] where powerful tools such as localisation [99] are available. Another 
approach studied extensively is based on the AdS/CFT correspondence [83, 100-107] where equilibration is thermal. 
An interesting open question is: if a thermal ensemble is described in this context by a black hole, what is the gravity- 
dual of a GGE ensemble? It is conjectured [67, 108, 109] that it corresponds to a higher spin black hole, whereas in 
classical gravity black holes can only have a small number of additional conserved charges (mass, angular momentum 
and electric charge) in contrast to the infinite number associated with integrability. It would be interesting to use our 


15 


approach in order to study such questions and compare with the other available methods. 
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Appendix A: Preliminary definitions 


We consider the Klein-Gordon field described by the hamiltonian (1). The time-evolved field operator <j> in the 
Heisenberg picture is 






Ak-i 


L d /2 


and similarly the conjugate momentum 7r is 


x/2 E, 


: (nfce 


-iE k t + 0 1 \_ k e+iE h t) 


(Al) 


< X A) = ^2 eik '^k(t) 


H) 

L d ! 2 


* \l^(a k e 


0 -iE k t _ t +iE k t\ 


(A2) 


Here L is the system size and we assume periodic boundary conditions, so that the momenta are given by k = 2nn/L 
with n integer. 

The creation-annihilation operators are 


/ Eb ~ % 

a k = \J -y^(O) + -j==TT k { 0), 


a -k = Y - ~^fk(Q)- 


The field (j>(x',t) automatically satisfies the equations of motion, which in momentum space read 

4>k{t) + El<f> k (t) = 0, 


with the obvious solution 


4 >k{t) = 4>k(0) cos Et~t + 7tfc(0) 


sin Ekt 
Ek 


(A3) 


(A4) 


(A5) 


Before we proceed further, we also report an identity that will be useful later in expressing the time averaged values 
of field correlations in terms of their initial values in coordinate space. From (A3) the operator n(k) + n(—k) can be 
written in terms of the field <f> and its time derivative (f> = tt as 


n(k) + n(—k) = a' k ak + aL fc a_ fc 
Ek^k^—k T Ek'E—k/Ek 1, 


(A6) 


where in the last step we used the commutation relations [</>*,, 7r q ] = i5k- q and = [^k^q] = 0. 
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Appendix B: The one-point correlation function 

Using the above expressions for the evolved fields (Appendix A), we can easily find the evolution of the one-point 
function. This is only non-zero if its initial value <po = {<P{ x i 0)) is non-zero (as for example in the case of the ground 
state of a QFT with broken symmetry). Indeed, from (A5) we obtain 

CfV) = {<t>(x;t)) = -^V fe ' x (0 fc (O))cos E k t 

k 

= </>o cos rat, (Bl) 

since 

(MO)) = J^J Ld d d x e~ ik ' x {4>{x] 0)) = 00 i d/2 4 d o, (B2) 

and (7Tfc(0)) = (d/dt)((j)k(t))\t = o = 0 as a consequence of the continuity of the time derivative of the field at t = 0 and 
the fact that before the quench there was no time evolution, since the initial state was an eigenstate of the pre-quench 
hamiltonian. From (Bl) we see that the evolution of the one-point function is trivially oscillating between the values 
±00- 

In the above and in what follows we make the assumption that the initial state is translationally invariant. Without 
loss of generality, from now on we will also assume that the initial value of the field is zero, i.e. 0o = 0. In the opposite 
case, we can always redefine the field 0 to be the difference from the above one-point function, therefore focusing on 
the quantum fluctuations about this value. 


Appendix C: The two-point correlation function 


Next we calculate the two-point function at equal times in terms of expectation values in the initial state 

C[ 2 \x,y) = ( 0 (aqf) 0 (y;f)) 

_ _L \ 1 _ ^ p ik 1 x+ik 2 y 

~ Ld tk 2 -J 2^7 

x {a kl a k2 )o e - i{E ^ +E ^ )t + (aL fci a fc ,) 0 e +z{E ^~ E ^ )t 

+ (afcial fc2 )o e-iCBur-E**)* + (aL fci aL fc2 ) 0 e +i{E ^ +E ^ ]t 

where the index zero means that the expectation values are calculated in the initial state. Using our assumption that 
the initial state is translationally invariant, the above expectation values are nonzero only when the momenta k i, k- 2 
are equal or opposite. We therefore have 


C[ 2 \x,y) = 0 V -i- e lk - {x ~ y) 

t V ’ L fc 2E k 


X 


1 + F(-k) + F(k) + G(-k) e 


— 2 iEkt 


+ G*(—k) e 


+2 iEkt 


(Cl) 


where G(k) and F(k) are functions that depend on the particular initial state: F{k ) expresses the momentum 
distribution in the initial state 


F ( k ) = (alafc) o = (n(k)) o, (C2) 

and it is a real function, while G(fe) gives the probability amplitude for the presence of a pair of opposite momentum 
particles in the initial state 


G(k) = {a- k a k ) 0 , (C3) 

and it is an even function. We will also assume that the initial state is rotationally invariant, in which case the above 
are functions of the norm k = \k\ but not of the direction of the momenta k. 
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In the thermodynamic limit L —> oo, the sum in the above expression becomes an integral. Using the translational 
and rotational symmetry of the initial state we can further simplify the above expressions, since the functions F and 
G are function of | fe | only. We now take the large time limit t —> oo of the correlation function. In the case of massive 
post-quench dispersion relation, Ek = Vk 2 + m 2 with m ^ 0, the stationary phase method shows that the oscillating 
time dependent terms vanish and therefore the two-point function equilibrates. In the massless case the large time 
asymptotics of the above correlation function is sensitive to the behaviour of the function G{k). Typically, if G(k ) 
decays sufficiently fast at large k and does not have a pole at k = 0 but possibly other singularities elsewhere in the 
complex fc-plane, then the time dependent terms decay for large times and the correlation function equilibrates. 

As long as the two-point function equilibrates, its stationary value is equal to its long time average value which is 

° {2) (r) = ^Y / ^-e ikr (l + 2F(k)). (C4) 


Note that it depends solely on the initial expectation values of the occupation numbers of the momentum modes 
F{k) = (n(fc)) 0 . 

Using the identity (A6) in (C4) we can rewrite it as 

C (2) (r) = X - (E< 2) (r) + J d d s H(s)C < 2) (r - a)) , (C5) 

L d 


“ ( 2 ) 

where the functions FI and Cq have been defined in (23) and (24) 

d d k e ik x 
(2 n) d E 2 ■ 


H(x) = J 


Co\x) =(7^(0,*), 


( 2 )/ 


and 


C^\x) = (7r(0,0)7r(a;,0)) 


d_d_ 
dti dt 2 


<0(O,ii)0(®,i 2 )) 


1 1 =£2=0 


(C6) 

(C7) 


(C8) 


Appendix D: The four-point correlation function 


We proceed to the calculation of the four-point function. The equal time four-point function can be expressed as 

C < f i {x 1 ,x 2 ,xz,x i ) = (cj)(xi;t)(f>(x 2 ; t)4>(x 3 -, t)(f>(x<±-, t)) 

= J_ V 1 n iT,U ki *i 

T2d Z_^ 


ki,k2,k3,k4 4 y/nU E ki 


E n< 

a11 { CT i} \ i—\ 

<Ti=± 


Gi) 

—(Tiki 


p i Ylt= 1 


(Dl) 


where we used the compact notation a = a^, and a^ ^ = a^. 

As before, assuming that the initial expectation values are suitably smooth functions so that the oscillating terms 
decay for large times, the four-point function equilibrates and its stationary value is given by the terms that are 
not accompanied by time-oscillating phases, i.e. those satisfying the condition ViEki = 0 with cq = ±. If the 

four-point function does not equilibrate, then this expression gives the long time average. Moreover, by translational 
invariance of the initial state, the above expectation values are non zero only if their total momentum is zero ki = 
0. Taking into account that the momenta k l are quantised, i.e. each of their components is an integer multiple of 
27r/L and that the first equation is not a linear equation (since Ef, = \J\k\ 2 + m 2 ), the terms that give non-vanishing 
contributions in the thermodynamic limit are those with equal number of a and a 1 operators, i.e. those for which 
XEi °i— 0 (which are 4!/(2!2!) = 6 in number, out of 2 4 = 16) and for which the momenta of the a operators match 
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with those of the a 1 operators in pairs of opposite values. All other solutions in the set of integer numbers are sparse 
and do not contribute in the thermodynamic limit. 

Considering all possible permutations we finally obtain an expression in terms of the initial expectation values of 
the occupation numbers n(k) and their products 


C (4) (x 1,X 2 ,X 3 ,X 4 ) 


1 1 
2L2d 


E 


S(k,p; x 4 ,x 2 , x 3 , x 4 ) 

4 J^k^p 


X 


[n(k)n{p)) o + (n(fe)) 0 


where the function S(k,p;x 1 , 2 : 2 , * 3 , * 4 ) is defined in (29) 


S(k,p-x 1 ,x 2 ,x 3 ,x 4 )= e ifc ^ 2 - Xl)+lp(x4 - X3) . 

all perm.s 
of 1,2,3,4 


(D2) 


(D3) 


In the thermodynamic limit we replace L~ 2d ► (2ir)~ 2d f d d kd d p. 

Finally, we write the last expression in terms of coordinate space initial correlation functions. The function 
S(k,p\xi 1 X 2 ,x 3 l x A ) is even in both k and p and symmetric under their interchange. Therefore we can replace 
the expression ((n(fc)n(p ))0 + (n(k )) 0 + |) in the sum by 

\ (0( fe ) + n (~ k ) + !) (n(p) + n(-p) + 1)) 0 = 

= ^((n{k)n(p)) 0 + {n(—k)n(p )) 0 + (n(k)n(-p)) 0 
+ (n(-k)n{-p)) 0 + (n(k )) 0 + (n(-k )) 0 + (n(p )) 0 = 

+ (n(-p)) 0 + l), (D4) 

and use the identity (A6) to write (D2) as 


1 1 


C (A \x i,x 2 ,x 3 ,x A ) = — -^yS(k,p-,x 1 ,x 2 ,x 3 ,x 4 ) 


32 L 2d 
1 /. 


k,p 


X ^ (^kty—k^p^—pj^ “I - -^2 \^k'^—k4 ) p4 ) —Pi ^ 

+ 2 k'kp'k -+ E 2 E 2 0 


(D5) 


Notice that the large time value of the four-point correlation function depends solely on four-point initial correlations, 
not on lower order correlations of the initial state. 

The last step is to write the correlation functions in coordinate space, thus obtaining 


C^\x ll x 2 ,x 3 ,x 4 ) = 


- y 

39 


32 ^ J J L 2d 

all perm.s , 
of 1,2,3,4 L d 


d d sd d r 


Cq 4) (s, s + xi - x 2 , r, r + x 3 - x 4 ) 


+ J d d s' H(s')Cq 4 \s, s' + s + x 4 - x 2 , r,r + x 3 - x 4 ) 

L d 

+ J d d r H(r')(C ( 0 4) (s, s + x 4 - x 2l r , r’ + r + x 3 - x 4 ))” 


+ 11 d d s'd d r' H(s l )H(r') C^\s, s' + s + x 4 — x 2 ,r 1 r' + r + x 3 — x 4 ) 


L d 


(D6) 
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where and C g^ are defined in (31) and (32) 

Cq A) {x ll x 2 ,x 3 ,x 4 ) = (n(x 1 )Tr(x 2 )^(x 3 )<l)(x 4 )) 0 (D7) 

and 

C < o\x 1 ,x 2 ,x 3 ,x 4 ) = (7r(x 1 )Tr(x 2 )Tr(x 3 )TT(x 4 )) 0 (D8) 

and the function H(x) has been defined in (23). 


Appendix E: Spatial averages of initial correlation functions 


We will show, using the cluster expansion and the cluster decomposition principle, that the spatial weighted averages 
of initial correlation functions entering in (30) are determined by the two-point functions only, provided that the weight 
function H(r ) decays exponentially with the distance r, or is, more generally, an integrable function in the whole 
infinite space. This argument is based on a multi-dimensional generalisation of the trivial fact that the infinite average 
value / of any function f(x) equals its value at infinity, provided that it is finite, i.e. 

f= \imL~ 1 [ dxf(x) = /(+oo). (El) 

L-yoo J 0 

Let us focus on the first term in the integral of (30), which is essentially a spatial average of the initial four-point 
correlation function Cq^ with respect to two coordinate variables s and r. This is then equal to its value at |s| —> oo 
and |r| —► oo, i.e. when the four points are separated in two pairs infinitely far from each other. According to the 
cluster decomposition principle, in this limit the correlation function tends to the disconnected expression, i.e. the 
product of two two-point correlation functions 


lim 

L—> oo 



L d 


d d sd d r (4) 
L id °o 


(s, S + Xi 


= C 3 \xi - x 2 )Cq 2 \x 3 - x 4 ). 


x 2 ,r,r + x 3 - x 4 ) 


(E2) 


In the above we used the fact that the initial expectation value of the field </> vanishes, {<j>(x))o = do = 0. Similar 
results hold for all other terms in the integral of (30). For the last term, for example, we find 


lim J J ^ r JJ d d s'd d r' H{s')H(r') ((^(s, s' + s + x 4 - x 2 ,r,r' + r + x 3 - x 4 ) 

L d L d 

= jj d d s'd d r' H(s')H(r') Cq 2 \x 1 - x 2 + s')Cq 2 \x 3 - x 4 + r'). (E3) 


Note that the last expression is always finite for a massive post-quench hamiltonian, because the function H(r) decays 
exponentially with the distance r. 

More explicitly, the general cluster expansion for the C'g 4 '* correlation function is 

Co 4 ’ (Sl, S 2 , S 3 , S 4 ) = C™ (Sl, S 2 ) & 2) (S 3 , Si) 

+ Co 2) (si,S3)Co 2) (s 2 ,s 4 ) + C'q 2 ^ (si, s 4 ) Cq 2 ' 1 (s 2 ,s 3 ) 

+ ( s l> s 2, s 3> s 4:) (E4) 

where Cg 4 ^ is the connected initial four-point function. Applying this formula and taking into account the translational 
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invariance of the initial state, we can write the integrand of the first term in (30) as 

Cq 4) (s, s + x 4 - x 2l r, r + x 3 - x 4 ) 

= Cq 2 \xi - x 2 )C^\x 3 - x 4 ) 

s v— 

type I 

+ Cq 2 \v - s)Cy 2 \r - s + x 3 - x 4 + x 2 - Xi) 

"---' 

type II 

+ C ( 0 2 \r - s + x 3 - x 4 )C^ 2 \r - s + x 2 - x 4 ) 

S -v- y 

type II 

+ Cq 4 c '(0, x 3 - x 2 ,r - s,r - s + x 3 — x 4 ) . (E5) 

S -v- ^ 

type III 

As before we omitted terms involving (0(ai)) 0 = 4 >o since the latter is assumed to be zero, which therefore means that 

( 2 ) ( 2 ) 

the connected two-point function is C a = 6q . 

We distinguish, as indicated in (E5), three types of terms: 

( 2 ) 

• type I is a product of two-point functions Cq ' whose arguments do not depend on the integration variables s 
and r, 

( 2 ) 

• type II are all other terms that are products of two-point functions Cq ' too, but do depend on the integration 
variables, 

• type III is the connected four-point function Cq 4 ^. 

The cluster decomposition principle requires that the connected correlation functions Cq 2 ^ = Cq 2 ' 1 and Cq 4 ^ decay 
to zero when the distance between their points increases. We then readily see that, substituting (E5) into (30) and 
performing the integration over s and r, only the type I term gives a contribution of order L 2d as required in order 
to survive in the thermodynamic limit, while all the rest give contributions of lower order in L. 

Similarly, for the C correlation function we have 

Cq 4) (s, s 1 + s + xi - x 2 , r,r' + r + x 3 - x 4 ) = 

= C^ 2) (aii - * 2 + s') C^ 2) (x 3 - x 4 + r') 

"-v-' 

type I 

+ Cq 2 ' 1 (r - s ) Cq 2) (r 1 — s' + r — s + x 3 - x 4 + x 2 - x{) 

s -v-' 

type II 

+ Cq 2 ' 1 (r' + r - s + a: 3 — x 4 ) (r - s - s' + x 2 - aq) 

s -v-' 

type II 

+ C’o 4 c (0, aci - ai 2 + s', r - s, r' + r - s + x 3 - x 4 ), (E6) 

s ---V-' 

type III 


and analogously for Cq 4 '. In the above, the type I term depends only on s' and r' but not on s and r as the other 
terms. 

If the post-quench hamiltonian is massive, then similar conclusions to the ones after (E5) hold for the contribution 
of the above terms in the thermodynamic limit after all integrations have been performed. Indeed, in this case the 
function H(r ), given by (23), decays exponentially with the distance, therefore once again only the type I term gives 
a contribution that scales as L 2d . 
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Summing up all surviving terms, (30) reduces to 


C (4) (x 1 ,x 2 ,x 3 ,x 4 ) = ^ ^2 

all perm.s 
of 1,2,3,4 

^Cq 2 \xi - x 2 ) + j d d s'H(s')CQ 2 \x 1 - x 2 + s')^j 

x ^C^ 2) (x 3 - x 4 ) + J d d r'H (r , )C l Q 2 ' ) ( x 3 - x 4 + r')^ , (E7) 

which, using (22) for the two-point correlation function, can be written as 

C ( - 4 \x ll x 2 ,x 3 ,x 4 ) 

= i V] (7 (2) (a;i - x 2 )C { - 2 \x 3 - x 4 ) 

o — 

all perm.s 
of 1,2,3,4 

= C (2 \x i - a: 2 )C ,(2) (a:3 - * 4 ) + [2 *+ 3] + [2 4], (E8) 

This is identical to Wick’s theorem expansion for the four-point function. 

If, on the other hand, the post-quench hamiltonian is massless, then H(r ) does not decay exponentially and therefore 
it is no longer certain that the asymptotic expressions of the integrals scale with L as they should in order to give 
a finite (non-divergent) contribution in the thermodynamic limit. More specifically, if H(r) diverges with L (as 
in d = 2) then the same is certainly true for the long time averaged correlation functions, which therefore do not 
equilibrate. Even if H{r) does not diverge with L but decays algebraically with the distance as r~ 2 or slower, then it 
is not an integrable weight function and the long time averaged functions may diverge with L, depending on the large 
distance behaviour of the initial correlation functions. However it is still true that the dominant contribution in the 
thermodynamic limit is given only by the type I terms above, since the cluster decomposition principle ensures that 
the spatial averages of all other types of terms scale slower with L. Therefore (E7) above is still applicable at leading 

order in L. We denote this by replacing the equality symbol by meaning that only the leading asymptotic 

expression is kept, while the omitted subdominant terms may also diverge with L but slower. Possible exceptions to 
this conclusion may arise only in cases where the initial correlation functions decay algebraically with exponents that 
have special values, such that type III terms do contribute to the leading-order in the thermodynamic limit. 


Appendix F: The GGE predictions for the correlation functions 


We will calculate the GGE prediction for the two-point correlation function. The GGE is given, always in the 
thermodynamic limit, by the (non-normalised) density matrix 

PGGE = exp J d d k /3 (fc)n(fc)^ , (FI) 

where the Lagrange multipliers /3(k) are defined through the requirement that the values of the charges n(k) in the 
GGE are equal to their initial values 


(n(fe)) gge = (n(fc)) 0 . 

The GGE value of the two-point function is then 

C'gge (*.y) = (H x )Hv)) gge 

__ e ik ( x -y) ^(oL fc a_fc)GGE + ( a fc a l) G GE^ 


d d k 

(2^p2 E k 


j Wf^E/ k ' (X ~ V) (1 + < n (" fc » GGE + H fc )>GGE) 


(F2) 


(F3) 


Since the GGE, as expressed in (4), is diagonal in the momentum basis, it is a gaussian ensemble and therefore 
Wick’s theorem for the expectation values of multi-point correlation functions is valid in it. This means that the 
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prediction for the four-point function is the sum of all possible disconnected two-point function contractions, that is 

^GGE^l’ *2) x 3’ X/ 0 

= ®2)Cgge( :e 3i *4) + [2 o 3] + [2 o 4] (F4) 

= r d d kd d p i , k . (x ,- X9 . )+lv . (x ,- Xi) 

J ( 2 it)™ 4 E k E p 

X ((n(— fe))cGE + (n(fe))GGE + 1) 

X ((n(-p))cGE + (n(p )}GGE + 1) 

+ [2*4 3]+ [2*4 4], 
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